LOTKA_VOLTERRA

Overview

The LOTKA_VOLTERRA function numerically solves the Lotka–Volterra predator-prey equations, a foundational model in mathematical biology that describes the dynamics of two interacting species—one as predator and one as prey. Developed independently by Alfred J. Lotka (1925) and Vito Volterra (1926), this system of coupled first-order nonlinear ordinary differential equations remains widely used in ecology, economics, and systems dynamics.

The model is defined by the following equations:

\frac{dx}{dt} = \alpha x - \beta xy

\frac{dy}{dt} = \delta xy - \gamma y

where x represents the prey population density, y represents the predator population density, \alpha is the prey birth rate, \beta is the predation rate coefficient, \delta is the predator reproduction rate per prey consumed, and \gamma is the predator death rate. The model produces characteristic oscillatory dynamics: prey populations grow exponentially in the absence of predators, predators increase when prey is abundant, and both populations cycle through peaks and troughs over time.

This implementation uses SciPy’s solve_ivp function to perform the numerical integration. The function supports multiple integration methods including explicit Runge-Kutta schemes (RK45, RK23, DOP853) for non-stiff problems and implicit methods (Radau, BDF, LSODA) for stiff systems. The default RK45 method uses a fifth-order Runge-Kutta formula with fourth-order error control. For more details, see the SciPy integration documentation and the SciPy GitHub repository.

The classic predator-prey oscillations predicted by this model have been observed in natural populations, most notably in the historical Hudson’s Bay Company fur trading records showing coupled fluctuations between lynx and snowshoe hare populations. For a comprehensive mathematical treatment, see the Wikipedia article on Lotka–Volterra equations.

This example function is provided as-is without any representation of accuracy.

Excel Usage

=LOTKA_VOLTERRA(prey_initial, predator_initial, alpha, beta, delta, gamma, t_start, t_end, timesteps, lv_method)

prey_initial (float, required): Initial number of prey.
predator_initial (float, required): Initial number of predators.
alpha (float, required): Prey birth rate.
beta (float, required): Predation rate coefficient.
delta (float, required): Predator reproduction rate per prey eaten.
gamma (float, required): Predator death rate.
t_start (float, required): Start time of integration.
t_end (float, required): End time of integration.
timesteps (int, optional, default: 10): Number of timesteps to solve for.
lv_method (str, optional, default: “RK45”): Integration method.

Returns (list[list]): 2D list with header [t, Prey, Predator], or error message string.

Examples

Example 1: Demo case 1

Inputs:

prey_initial	predator_initial	alpha	beta	delta	gamma	t_start	t_end	timesteps
40	9	0.1	0.02	0.01	0.1	0	20	10

Excel formula:

=LOTKA_VOLTERRA(40, 9, 0.1, 0.02, 0.01, 0.1, 0, 20, 10)

Expected output:

t	Prey	Predator
0	40	9
2.222	28.93	15.7
4.444	15.9	20.57
6.667	7.737	21.2
8.889	3.914	19.21
11.11	2.212	16.43
13.33	1.416	13.68
15.56	1.018	11.25
17.78	0.808	9.19
20	0.6975	7.482

Example 2: Demo case 2

Inputs:

prey_initial	predator_initial	alpha	beta	delta	gamma	t_start	t_end	timesteps	lv_method
40	9	0.1	0.02	0.01	0.1	0	20	10	RK23

Excel formula:

=LOTKA_VOLTERRA(40, 9, 0.1, 0.02, 0.01, 0.1, 0, 20, 10, "RK23")

Expected output:

t	Prey	Predator
0	40	9
2.222	28.93	15.7
4.444	15.9	20.57
6.667	7.737	21.2
8.889	3.914	19.21
11.11	2.212	16.43
13.33	1.416	13.68
15.56	1.018	11.25
17.78	0.808	9.19
20	0.6975	7.482

Example 3: Demo case 3

Inputs:

prey_initial	predator_initial	alpha	beta	delta	gamma	t_start	t_end	timesteps
50	10	0.1	0.02	0.01	0.1	0	40	10

Excel formula:

=LOTKA_VOLTERRA(50, 10, 0.1, 0.02, 0.01, 0.1, 0, 40, 10)

Expected output:

t	Prey	Predator
0	50	10
4.444	13.67	26.57
8.889	2.229	22.36
13.33	0.6646	15.14
17.78	0.3457	9.911
22.22	0.2635	6.438
26.67	0.2583	4.176
31.11	0.298	2.71
35.56	0.3822	1.764
40	0.5246	1.154

Example 4: Demo case 4

Inputs:

prey_initial	predator_initial	alpha	beta	delta	gamma	t_start	t_end	timesteps	lv_method
30	5	0.15	0.025	0.015	0.12	0	10	10	RK45

Excel formula:

=LOTKA_VOLTERRA(30, 5, 0.15, 0.025, 0.015, 0.12, 0, 10, 10, "RK45")

Expected output:

t	Prey	Predator
0	30	5
1.111	29.96	7.231
2.222	27.8	10.27
3.333	23.51	13.82
4.444	18.04	17.12
5.556	12.81	19.35
6.667	8.706	20.22
7.778	5.874	19.95
8.889	4.038	18.94
10	2.872	17.55

Python Code

from scipy.integrate import solve_ivp
import numpy as np

def lotka_volterra(prey_initial, predator_initial, alpha, beta, delta, gamma, t_start, t_end, timesteps=10, lv_method='RK45'):
    """
    Numerically solves the Lotka-Volterra predator-prey system of ordinary differential equations.

    See: https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.solve_ivp.html

    This example function is provided as-is without any representation of accuracy.

    Args:
        prey_initial (float): Initial number of prey.
        predator_initial (float): Initial number of predators.
        alpha (float): Prey birth rate.
        beta (float): Predation rate coefficient.
        delta (float): Predator reproduction rate per prey eaten.
        gamma (float): Predator death rate.
        t_start (float): Start time of integration.
        t_end (float): End time of integration.
        timesteps (int, optional): Number of timesteps to solve for. Default is 10.
        lv_method (str, optional): Integration method. Valid options: RK45, RK23, DOP853, Radau, BDF, LSODA. Default is 'RK45'.

    Returns:
        list[list]: 2D list with header [t, Prey, Predator], or error message string.
    """
    # Validate input types
    try:
        prey0 = float(prey_initial)
        predator0 = float(predator_initial)
        a = float(alpha)
        b = float(beta)
        d = float(delta)
        g = float(gamma)
        t0 = float(t_start)
        t1 = float(t_end)
        ntp = int(timesteps)
    except Exception:
        return "Invalid input: all parameters must be numbers."
    if t1 <= t0:
        return "Invalid input: t_end must be greater than t_start."
    if ntp <= 0:
        return "Invalid input: timesteps must be a positive integer."
    if prey0 < 0 or predator0 < 0:
        return "Invalid input: initial populations must be non-negative."
    if a < 0 or b < 0 or d < 0 or g < 0:
        return "Invalid input: rate parameters must be non-negative."
    # Validate lv_method
    allowed_methods = ['RK45', 'RK23', 'DOP853', 'Radau', 'BDF', 'LSODA']
    if lv_method not in allowed_methods:
        return f"Invalid input: lv_method must be one of {allowed_methods}."
    # Create time array for evaluation
    t_eval = np.linspace(t0, t1, ntp)

    # Lotka-Volterra ODE system
    def lv_ode(t, y):
        prey, predator = y
        dprey_dt = a * prey - b * prey * predator
        dpredator_dt = d * prey * predator - g * predator
        return [dprey_dt, dpredator_dt]
    # Integrate
    try:
        sol = solve_ivp(
            lv_ode,
            [t0, t1],
            [prey0, predator0],
            method=lv_method,
            t_eval=t_eval,
            vectorized=False,
            rtol=1e-6,
            atol=1e-8
        )
    except Exception as e:
        return f"solve_ivp error: {e}"
    if not sol.success:
        return f"Integration failed: {sol.message}"
    # Prepare output: header row + data rows
    result = [["t", "Prey", "Predator"]]
    for i in range(len(sol.t)):
        t_val = float(sol.t[i])
        prey_val = float(sol.y[0][i])
        predator_val = float(sol.y[1][i])
        # Disallow nan/inf
        if any([
            t_val != t_val or t_val == float('inf') or t_val == float('-inf'),
            prey_val != prey_val or prey_val == float('inf') or prey_val == float('-inf'),
            predator_val != predator_val or predator_val == float('inf') or predator_val == float('-inf')
        ]):
            return "Invalid output: nan or inf encountered in results."
        result.append([t_val, prey_val, predator_val])
    return result

Overview

Excel Usage

Examples

Python Code

Online Calculator